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ABSTRACT 

Stochastic knapsack problem originally was a versatile model for controls in telecom- 
munication networks. Recently, it draws attentions of revenue management community by 
serving as a basic model for allocating resources over time. We develop approximation 
schemes for knapsack problems in this paper, a system of nonlinear but solvable partial 
differential equations and stochastic partial differential equation are shown to be the limit of 
the process that following the optimal solution of the stochastic knapsack problem. 
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1 Introduction 

Stochastic knapsack problem(a.k.a. dynamic knapsack problem) is a typical mathematical 
model for sequential resource allocation, see, e.g. (?). Similar model was proposed also for 
queueing control problems, see, e.g. (?), (?) and (?). In the 1990's, this type of models were 
intensively studied in the analysis and design of telecommunication networks, see, e.g. (?) 
and (?). Lately, the stochastic knapsack problem drew considerable attentions in the opera- 
tions research community due to their applications in revenue (yield) management, see, (?) 
and references therein. The solution of the problem can generally obtained through dynamic 
programming. However, these stochastic knapsack problems are usually only components 
of a larger optimal control problem. For example, in many revenue management problems, 
stochastic knapsack serves as a realistic model for pricing activities. In order to achieve over- 
all financial objectives, a pricing model have to be tied with manufacturing and supply chain 
management. In these cases, constantly checking the look-up table generated by dynamic 
programming will not be feasible. Efficient approximations in more explicit forms are hence 
needed. In this paper, we intend to study the first and second order approximations in concise 
form. 



The methodology we follow is basically the fluid and diffusion approximation approaches that 
derived from functional law of large numbers and functional central limit theorems. These 
methods have been applied very successfully in the study of queues and queueing networks, 
see, e.g. (?). Our study can be viewed as a further enrichment to this body of research. In 
queueing context, the limit process are usually functional of simple differential equations and 
(reflected) linear differential equations. The limit we obtain is a solution of Monge-Ampere 
equation, one of the most important nonlinear partial differential equations. Monge-Ampere 
equation was linked to crucial quantities in geometry, was extensively studied by researchers 
in various branches of mathematics. We also hope that this could lead to further interactions 
between related fields. 

The specific setup of the problem under consideration is the following, at time t = 0, we have 
W units of resource available, at any time t = 1, • • • ,T, a request(demand) arrives in the form 
of a bivariate random vector, (P t , Q t ), where the two components represent the unit offer price 
and the required quantity; at each time t, after the realization of the request, we decide whether 
to accept or reject it to maximize the average revenue collected by time T. Furthermore, we 
assume that, 

• Suppose (P t , Q t ), for each t (and i.i.d. across t) follows a joint (discrete) distribution: 

P[Pt=p i ,Qt = q j ]-=O ij , i = 0,l,...,£; j = l,...,k. (1.1) 

Here, p = 0, and 9 ■.= X^#oj represents the probability that there is no arrival (re- 
quest/demand) in a given period; 

• No partial fulfillment is allowed. 

The decision of accept/rejection at each time t, given that the available resource level is d, is 
determined by the following stochastic dynamic programming, 

V(t,d) = V(t + l,d)[9 + e(d)] 

+ E E d ij-™ a *{Pi ( lj + V( t + 1 , d -qj),V(t + 1 ,d)}, (1.2) 

iy^O j.Qj<d 

where 

W : =E E 9 *- ( 1 - 3 ) 

i -• j-qj>d 

Clearly, the first term on the right hand side of (|1 .2) corresponds to the case of either no arrival 
or the request size exceeds the available inventory; whereas each term under the summation 
compares the two actions: accept (i.e., supply) the request, or reject it. If we supply the request, 
then we earn the revenue pjqj, and proceed to the next period with qj units less in the available 
inventory. 



At the last period, we have 

v(T,d) = J2 E e vm> ( 1 - 4 ) 

*v^0 j:qj<d 

since clearly the best action is to supply any possible requests using the remaining inventory. 
The optimal decision is hence, to accept the demand (P t , Q t ), if it satisfies, 

PtQt + V(t + l,d-Q t )>V(t + l,d); (1.5) 

and reject it otherwise. 

This problem arises from many applications involving sequentially allocating limit amount of 
resources to multiple classes of demands, hence has been widely studied. The existing re- 
search can be divided into two groups, one focuses on the structural properties of the value 
function, usually through Markov decision process and modularity argument, see, e.g. (?), (?); 
the other focus on a variety of heuristics policies and their performance, much attentions was 
drawn some applications in the field of revenue management and dynamic pricing, which is just 
another mechanism to implement the accept/rejection decision, see, e.g. (?). In this paper, we 
examine the asymptotic behaviors of the value function. Our approach is similar to the "fluid 
limits" and "diffusion limits" methods in queueing theory, where time and space are properly 
scaled to induce asymptotic behaviors of the systems. To be more specific, we scale both time 
and space, which is the amount we allow to allocate, by a constant n, then let n go to infinity, 
and examine the asymptotic behaviors of the value function under this situation. Equivalently, 
we consider that following two processes, 

V n (t, d) = -Vint, nd), V n (t, d) = -=V(nt, nd), 
n 

In the following two sections, we will exam the asymptotic behavior of V as n — > oo, we show 
that it converges to the solution of a boundary value problem of a homogeneous Monge- 
Ampere equation, give the general solution to the equation, and solve it for some special 
case; then we in order to characterize the rate of the convergence and the asymptotic random 
behaviors, we identify the stochastic differential equation that the limit of V n (t,d) satisfies. 
When the knapsack is characterized by more than one characters, we have a multi-dimensional 
stochastic knapsack problem. It is discovered recently to be an important model in the study 
of inventory management. Unlike the one dimensional case, the dynamic programming for the 
multi-dimensional stochastic knapsack problem is exponential with respect to the size of the 
problem. Approximations can certainly play even more important roles. Our fluid and diffusion 
approximation are carried over to this important model. A system of differential equation needs 
to be solved to get the fluid limit. 

The reset of the paper will be organized as the following, in Sec. H we will establish the fluid 
limit through probabilistic arguments, and present the solution to the Monge-Ampere equation; 
in Sec. EH diffusion approximations are established for the case of unit demand; then, the 
model is extended to the case of multi-dimensional stochastic knapsack problem, for both fluid 
and diffusion approximation. 



2 Convergence Results 



To facilitate the analysis, we will take a de tour. Instead of the value function, we study the 
stochastic process whose mean will achieve the value function. Let us denote X t ,d(s),s = 
t, t + 1, • • ■ , T to be the reward collected at time s following the optimal policy while starting at 
time t with d units available. Apparently, V(t,d) = E[X ttd (T)]. Define, 

V n (t,d) = -X ntM {nT). (2.1) 
n 

We expect that a law of large number type result exists, so that V n (t,d) can converge to 
V(t,d) for each (t,d). Moreover, we hope to extract the dynamic it satisfies, thus give us the 
relationship we expect from V(t, d). 
To prove the convergence, we need, 

Lemma 1 . Let V n (t, d) be the random variable defined above, then, 

VaY[X ntM (nT)} = 0(n), as n -> oo. (2.2) 

Proof We prove this lemma by induction. First, when n = 1, observe that when d approach 
infinity, the optimal policy will be just accept any arrivals, then the variance will be constant. 
Therefore, we have a maximum for the variance for all d. Denote M to be twice that number, 
we want to show that for each n, 

\/ar[X ntnd (nT)] ww , 
n 

Suppose it holds up to n. Then for the case of n + 1, apparently, we need, 

Var[X (w+1)f|(w+1)d ((n + l)r)] < M ^ 
n + 1 

To see this, conditional upon the arrivals from (n + l)i + l to (n + l)t + (T-t), whose cr-algebra 
can be denoted as Ffa+yt+iT-t) we have, 

Var[X (n+1)t)(n+1)d ((n + 1)T)] = E[Var[Jf (n+1)ti(n+1)d ((n + l)T)|jr (n+1)t+(T _ t) ]] 

+Var[E[X (n+1)ti(n+1)d ((n + l)T)|^ (n+1)t+(T _ t) ]], 

in which, 

E[Var[X (n+1)ti(n+1)d ((n + l)T)\F {n+1)tHT _ t) }} < nM, 
from inductive assumption, and 

Var[E[X {n+1)ti(n+1)(i ((n + l)T)|^ {n+1)t+(T _ t) ]] < M. 
due to the definition of M. □ 

For the convergence proof, we adapt the classical proof of the strong law of large number, see, 
e.g. (?). Equation (|2.2) > and Markov inequality imply, 



n=l 



X u n t,u n d{u n T) - E[X Unt ^ Und (u n T)] 



> e 



< oo, 



where u n = [a n \ for some a > 1. From Borel-Cantelli lemma, we know that X Unt . and {u n T)/u n 
converges almost surely, and let us denote the limit as u(t,d) for each (t,d). Meanwhile, for 

any u n <k < u n+ i, we have, 

u n X Unt)Und {u n T) X k t,kd(kT) u n+ i X Un+lt ^ Un+ld (u n+ iT) 



U n +1 U n K U n U n+ i 

This leads to, 

1 u ^ ■ r x kt,kd(kT) X kt kd (kT) 

—u(t,a) < limint < limsup : < au{t,d), 

a k k k k 

since it holds for any a > 1, we can conclude that the convergence with probability one. 
Next, we relax the integral assumptions on the time and the knapsack size. For any (t,d) e 
M+, define, X nt)nd {nT) = Xy nt ^\nt\{nT). Then it is easy to see that the convergence can be 
extended. Furthermore, the converge is uniform on compact set, for related details, see. e.g. 
(?). Therefore we have, 

Theorem 2. V n (t,d) -> u(t,d) almost surely as n -> oo. and u(t,d) satisfies the Monge- 
Ampere equation with proper boundary value. 



du(x,y) , ( du(x,y) s _ n 

dx ' 9\ d y ) — u > <2 3) 

u(X, y) = h(y), u(x, 0) = 0, (x, y) £ [0, X] x [0, Y] 



where g(-) is the loss function of the random variable PQ, i.e. g(x) •— E[Q[P - x} + ]. 

Proof The only thing left to proof is that the limit will satisfies the equation in (|2.3|) . Since 
the convergence is uniform on compact set, we can further conclude that, if u is continuously 
differentiable with respect to t and dona compact set, which imply uniform continuity for both 
u and its derivatives, thenQ 

(nT) 

X ndM {nT) - X nd , nt+1 (nT) = a a ► ^> a - s - 

1/n ot 



du(t, d) 



X n d,nt(nT) - X nd+ i nt (nT) = ^ — ^ ► -J — .a.s. 

1/n ad 

This can be achieve from the following relationship. 

X ndM {nT) - X ndM+1 {nT) = [V{nd -Q,nt + l)+PQ- V(nd,nt + 1)]+ (2.4) 
This is the same as, 

Xnd,nt( nT ) Xnd,nt+i(nT) V (nd-Q ,nt+l) V(nd,nt+1) 

a j a = [PQ + 1 ] + - (2.5) 

n n 



'The following relations can be also derived through actions of Schwartz functions to define derivative in the 
distribution sense, hence, the solutions of the partial differential differential equation are weak solutions. However, 
in our case, the PDE has a classical solution, this direct approach is more intuitive. 



Since, both Xnd ^ nT ^> and V n (d,t), as its mean, converge to the same function uniformly on 
compact set, and according to regularity theorems on the first order partial differential equa- 
tions, both this function and its derivatives are uniformly continuous, see, e.g. (?), we can 
ensure the convergence of differences to the derivatives. (|2.5|) gives us the relation, 

du(x,y) du(x,y) 

— a — o ) = °- 

ox Oy 

We can draw the same conclusion from the boundary value. □ 

The equation (|2.3} is, of course, a nonlinear partial differential equation, which in general is 
extremely difficult to analyze. However, we can reduce this problem to a very well known 
equation. Just take derivative with respect to x, we have, 

UXX + g'(Uy)u xy = 0, 

similarly, take derivative with respect to y, we have 

U X y + g'(Uy)Uyy = 0, 

multiply the first one by u yy and the second one by u xy , then take difference, we have, u xx u yy - 
ul y = 0. This is a homogeneous Monge-Ampere equation. In general, Monge- Ampere equa- 
tion is an important mathematical subject in analysis and geometry. The general solution to 
the homogeneous Monge-Ampere solutions was obtained independently in (?) and (?). We 
will follow the basic setup and structure of the latter. In Fairlie and Leznov(1 995) (?), a general 
construction to the second class is obtained. Given any two arbitrary C 1 function R(-) and /(•), 
a solution can be uniquely determined. Here, the general procedure of the construction gives 
us the basic forms that the solution will take, g(-) and h(-) then will help us to determine the 
solution completely. We will also reveal the relations between the two sets of functions. The 
general solution to the homogeneous Monge-Ampere equation bears the following form, 

du dR du dR 

dx~ = { ^~^ , ~oy = ~dj' 

V- a* = /(£)• (2.6) 

For any continuous function R(£) and /(£). For our problem, we know that £+g(R')+g' (R 1 ) = 0, 
u(X,y) = h(y), these conditions can uniquely determine the two functions R(£) and /(£), 
hence, uniquely determine u(x, y). For example, when the batch distribution takes the following 
special form, g(x) = e~ Xx , = - 1 ° s( ^ l) , we know that 



u(x,y) 



l' V dR J 

' + ' ~dj 



!J(i 



in conjunction with the relation, 

d£ _ £ d£_ 1 

^ " x-r(zy dy~ - f(o-x- 

we can obtain u(x,y) from routine differential equation procedures for first order partial differ- 
ential equation, see e.g. (?). 



3 Second Order Approximations 



While fluid limits capture the dynamics of the system evolution, asymptotic characterization of 
the randomness of the system is also desired. Diffusion approximations, results of various cen- 
tral limit theorems and strong approximation theorems usually serve this purpose. In addition, 
diffusion approximations can also provide an estimation of the order of the system converging 
to the fluid limits. The process understudy differs from most of those appears in the queueing 
context is that the main process is a random walk with time and state dependency, furthermore, 
it depends upon the increments of another random walk. Standard functional central limit the- 
orem can not be applied directly. How ever in the case of unit demand, we can overcome this 
barrier by making use of results derived in (?). In (?), a large class of switching system are 
studied. Our system turns out to be one that "switching" at every time period. To apply the 
results in (?), it is require that the evolution function of the mean and variance are smooth, in 
fact Lipschitz should be enough. Those functions for our case are functions of the solutions 
to the Monge-Ampere equation, whose smoothness has been well established. Then, results 
in (?) will give us the stochastic differential equation(SDE) the limiting process will satisfy. 
However, the existence and uniqueness of the solution to the SDE were not established in 
(?). Meanwhile, in (?), the existence and uniqueness of the solution to this type of SDE with 
more general assumptions on the coefficients have been obtained. Hence, the SDE uniquely 
determine a stochastic process. 

To capture the randomness of X njd (T), we define yd,r(s) to be the number of units that are 
supplied to the demands up to time s. From the set up, we know that y d ,T{s) is a Markov chain. 
Moreover, we have, 

y d ,r( s + 1) = Vd,T(a) + X ({V(s + l,d- y d , T (s) - 1) +p - V(s + 1, d - y d , T (s)) > 0}), 

Z d:T (s) - Z d:T (s) = Pi[y d ,T{s) - Vd,T(s)}. 

In (?), random walks of this type have been studied. More specifically, Let {£„(«), a e M + } be 
independent random variables with al as parameters, and S n+ i = S n + £(S n ). The following 
theorem is a summary of results in (?), 

Theorem 3. ((?)) If there exists a function s(t) such that, 

\-S n (nt) - s(t)\ -► 0,a.s. 
n 

and £> n -»■ D(a), then, D n = Var[£ n (a)]/n then, ^ n {t) ■- n~ 1/2 [S n (nt) - ns(t)] weekly converge 
in D to the following stochastic differential equation, 

d 7 (t) = D(s(t)) l l 2 dw u 

provided that its solution exists and is unique. 



Following the similar functional strong law of large numbers in the previous section, we can 
see that s(t) satisfies, 

By the Lipschitz continuity of the solution to the Monge-Ampere equation, and existence and 
stability theorems of ordinary differential equation, we can see that s(t) can be uniquely de- 
termined and is Lipschtz continuous. Therefore, we can conclude that n _1 / 2 [F n (ni) - ns{t)} 
converge to a stochastic process that satisfies the following stochastic differential equation, 

dY t = V9(Mt, d - s(t)))[l - g(u d (t, d - s(t)))]dW t , (3.1) 

where dW t denotes a wiener integral with respect to a standard Brownian motion, and the 
revenue collected can be represented by 

dZ s = PdY s . (3.2) 

Remark: For the problems with general demand batch size, we conjecture that the process 
will weakly converge to a process that satisfies the following stochastic differential equation, 

dY t = y/F p>d (u d (t, d - Y t ))[l - F p>d (u d (t, d - Y t ))]dW t . 

Where F P:d is the distribution function for both price and the quantity, and is not necessary 
continuous. Hence, the limiting process is not a diffusion process any more. However, results 
in (?) guarantee the existence and uniqueness of the limiting process. This convergence will 
be the subject of future research. 



4 Multi-dimensional Stochastic Knapsack Problems 

In this section, we intend to deal with a multi-dimensional version of the stochastic knapsack 
problem. To be more specific, at time t = 0, we have a vector of m different types of prod- 
ucts available, {W l ,W 2 ,--- ,W m ), at any time t = 1,2, ■ ■ • ,T, the demand arrives now is 
characterized by a m + 1 tube, (P t , Q\, Q 2 , ■ ■ ■ , Q™), where P t describes the rewards and a 
m-dimensional characterization vector (Q\, Qf, ■■■ , Q™) describes the combination of the the 
m types products it requires. They are i.i.d random variables, with discrete distribution, 

P[P t = Pi, Ql =qj 1 ,Q 2 = q%,--- ,QT = CJ = OiJiJ*,- J m - 

Again, at each time t, we need to determine whether to accept or to reject the demand arrival 
to achieve the maximum average revenue at time T. Let us denote V(t, d 1 , d 2 , ■ ■ ■ , d m ) be the 
value function for the dynamic programming. Similarly to the one dimensional case, it satisfies 



the following recursive equation, 



V{t,d\d 2 ,--- ,d m ) 
= V(t + l,d\d 2 ,--- ,d m )[e + Q{d l ,d 2 ,--- ,d m )} 

■ mayi {Pi + V(t+ l,d l - qi,d 2 -q 2 ,--- ,d m - q m ), 

V(t + l,d\d 2 ,--- ,d m )}. (4.1) 

Here, 

@(d\d 2 ,--- ,dr) = P[uf =1 {Q k >d k }]. 

It is obvious that the dynamic programming is not so easy-to-solve as those in one-dimensional 
case. In fact, the computational efforts needed grow exponentially with respect to W = 
maxfc W k . Therefore, approximation will be much more desirable. 

Apply the sam scaling and argument, we can see that the similar convergence results hold, 

Theorem 4. V n (t, d 1 , d 2 , ■ ■ ■ ,d m ) -► u(t,d\d 2 ,- ■ ■ ,d m ) as n -► oo. and u(t,d\d 2 ,- ■ ■ ,d m ) 
satisfies the following boundary value problem, 

( du(xp,xi,X2,--- ,x n -i) , rit du du du \ _ n 

J 9x ^A&ci ' 9x2 ' " ' ' dx m -! I — u ' 

I u{xq, xi,x 2 ,--- , x m -i) = h(x 1 ,x 2 , ■■■ , x n -i),u(x ,0) = 0, 

where 

G(z±,z 2 , • • • , Zm-i) ■= E[P — [z\ + z 2 -\ h z m _i)] + . 

Hence, u satisfies the homogeneous Monge-Ampere equation Det(D 2 u) = 0. 

Proof: The only thing left is to show that u satisfy the Monge-Ampere equation. To see this, 
take derivative with respect to xq,x\, ■ ■ ■ , x m respectively on equation, 

du(x ,xi,x 2 , ■ ■ ■ ,x n -i) | du du du Q 

Oxq dx\ ' dx 2 ' ' dx m ^i 

we have, for k = 0, 1, • • ■ ,m 



du 2 ^ ^ du du du du 2 ^ ^ 

dx Q x k ^ 1 dx\ ' dx 2 ' ' dx m -i dxiXu 

regard this as a linear system for variables (1, Gi,G 2 ,-- , G m ), then in order for (|4~2]> to hold, 
the coefficient matrix, (D 2 u), must be singular. □ 

Again from (?), we have the general solution for n dimensional Monge-Ampere equation. Un- 
derstandably, it depends upon the selection of two m - 1 variate functions. Given arbitrary 
function - ,Cn-i), • • • ,Sn-i), 

u = R — ^j-JSj , Uj = iij , 



Q = (D 2 Ry 1 DL. 

Where DL and D 2 R denote the gradient vector and Hessian matrix of L and R respectively. 
Function G(-) and the boundary condition h(-) will then determine R and L. 
Define y^ T {s), k = 1, 2, ■ ■ ■ , m to be the number of units of characterizations that are supplied 
to the demands up to time s. There exists a function (s 1 (t) , s 2 (t) , ■ ■ ■ ,s m (t)) which is the 
solution to the following differential equation system, 

= G{u i (t,d 1 -s 1 ^)^ 2 -s 2 {t),--- ,d m -s m (t)), k = 1,2,- ■■ ,m. 

The regularity of the solution will then help us conclude that, 

fY£(nt) -ns k (t)\ 
V n^ 2 ) ' 

converges to a stochastic process determined by the following stochastic differential equation, 

dY t k = ^G k (l-G k )dW t k , (4.3) 
where, G k = G^^d 1 - s x {t),d? - s 2 (t), ■ ■ ■ ,d m -s m (t)). 



